// CSc 522, hw1
// multi-thread version of red/black grid computation
//Jinyan Huan (jguan1@email.arizona.edu) & Fengqiong Huang (huangfq@emailo.arizona.edu)

#include <math.h>
#include <stdio.h>
#include <sys/types.h>
#include <stdlib.h>
#include <pthread.h>
#include <sys/time.h>

// Global variables 
volatile int *arrive;
double **grid;
double *max_diffs;
int height;
struct timeval timer;
double start_time;
double end_time;
int num_threads, grid_size, num_iters;


// helper functions
void print_grid(double **mat) {
	int i, j;

	printf("The %d * %d grid is\n", grid_size, grid_size);
	for (i = 0; i <= grid_size + 1; i++) {
		for (j = 0; j <= grid_size + 1; j++) {
			printf("%lf ", mat[i][j]);
		}
		printf("\n");
	}
}

// Allocate grid
double ** allocate_grid(int size) {
	int i;
	double *vals, **temp;

	// allocate values
	vals = (double *) malloc(size * size * sizeof(double));

	// allocate vector of pointers
	temp = (double **) malloc(size * sizeof(double*));

	for (i = 0; i < size; i++) {
		temp[i] = &(vals[i * size]);
	}

	return temp;
}

// Barrier for id's thread using desemination barrier from class
void barrier(int id) {
	int j, look_at;
	for (j = 1;  j <= ceil(log(num_threads)); j++) {
		while (arrive[id] != 0)
            ;
		arrive[id] = j;
		look_at = (id + (int)pow(2, j-1)) % num_threads;
		while (arrive[look_at] != j) 
            ;
		arrive[look_at] = 0;
	}
}

// Perform red-black grid-computation for i'th thread
void compute(int id) {
	int first_row = id * height + 1;
	int last_row = first_row + height - 1;
	int j_start;
	int iter;
	int i, j;
	double my_diff;
	double temp, old_val;

	for (iter = 1; iter <= num_iters - 1; iter++) {
		// update red points
		for (i = first_row; i <= last_row; i++) {
			if (i % 2 == 1)
				j_start = 1;
			else
				j_start = 2;
			for (j = j_start; j <= grid_size; j = j + 2) {
				grid[i][j] = (grid[i - 1][j] + grid[i][j - 1] + grid[i + 1][j]
						+ grid[i][j + 1]) * 0.25;
			}
		}

        // synchronize the threads before updating the black points to make sure 
        //all red points have been computed
		barrier(id);

		// update black points
		for (i = first_row; i <= last_row; i++) {
			if (i % 2 == 1)
				j_start = 2;
			else
				j_start = 1;
			for (j = j_start; j <= grid_size; j = j + 2) {
				grid[i][j] = (grid[i - 1][j] + grid[i][j - 1] + grid[i + 1][j]
						+ grid[i][j + 1]) * 0.25;
			}
		}

        // synchronize the threads for the next round's update
		barrier(id);
	}

	// One more updates and keep track of max_diff
	// update red points

	my_diff = 0.0;
	for (i = first_row; i <= last_row; i++) {
		if (i % 2 == 1)
			j_start = 1;
		else
			j_start = 2;
		for (j = j_start; j <= grid_size; j = j + 2) {
			old_val = grid[i][j];
			grid[i][j] = (grid[i - 1][j] + grid[i][j - 1] + grid[i + 1][j]
					+ grid[i][j + 1]) * 0.25;
			temp = old_val - grid[i][j];
			if (temp < 0)
				temp = -temp;
			if (temp > my_diff)
				my_diff = temp;
		}
	}

    // synchronize the threads before updating the black points to make sure 
    //all red points have been computed
	barrier(id);

	// update black points
	for (i = first_row; i <= last_row; i++) {
		if (i % 2 == 1)
			j_start = 2;
		else
			j_start = 1;
		for (j = j_start; j <= grid_size; j = j + 2) {
			old_val = grid[i][j];
			grid[i][j] = (grid[i - 1][j] + grid[i][j - 1] + grid[i + 1][j]
					+ grid[i][j + 1]) * 0.25;
			temp = old_val - grid[i][j];
			if (temp < 0)
				temp = -temp;
			if (temp > my_diff)
				my_diff = temp;
		}
	}

	max_diffs[id] = my_diff;
    
    //no need to add barrier here since the join of all threads works like a barrier
}

// worker
void *worker(void *arg) {
	int id = *((int *) arg);
	compute(id);
}

int main(int argc, char** argv) {
	pthread_t *threads;
	double max_diff;
	int *p;
	int i, j;

	if (argc != 4) {
		printf("Usuage: %s <grid-size> <num-iterations> <num-threads> \n",
				argv[0]);
		exit(1);
	}
	grid_size = atoi(argv[1]);
	num_iters = atoi(argv[2]);
	num_threads = atoi(argv[3]);

	//validate program arguments
	if (grid_size <= 0 || grid_size % 8 != 0) {
		printf("Invalid grid_size: grid_size should be a multiple of 8\n");
		exit(1);
	}
	if (num_iters <= 0) {
		printf("Invalid num_iters: num_iters should be positive integer.\n");
		exit(1);
	}
	if (num_threads <= 0) {
		printf("Invalid num_threads: num_threads should be positive integer.\n");
		exit(1);
	}

	// Compute HEIGHT
	height = grid_size / num_threads;

	// allocate vector of pointers
	grid = (double **) allocate_grid(grid_size + 2);

	// Initialize the grid
	for (i = 1; i < grid_size + 1; i++) {
		for (j = 1; j < grid_size + 1; j++) {
			grid[i][j] = 0.0;
		}
	}
	// borders
	for (j = 0; j < grid_size + 2; j++) {
		grid[0][j] = 1.0;
		grid[grid_size + 1][j] = 1.0;
	}
	for (i = 0; i < grid_size + 2; i++) {
		grid[i][0] = 1.0;
		grid[i][grid_size + 1] = 1.0;
	}
	//print_grid(grid);

	// Allocate thread handles
	threads = (pthread_t *) malloc(num_threads * sizeof(pthread_t));

	// Initialize global arrive
	arrive = (int *) malloc(num_threads * sizeof(int));

    for(i = 0; i < num_threads; i++)
    {
        arrive[i] = 0.0;
    }

	// Initialize the max_diffs
	max_diffs = (double *) malloc(num_threads * sizeof(double));

	for (i = 0; i < num_threads; i++) {
		max_diffs[i] = 0.0;
	}

	// Start timer
	gettimeofday(&timer, NULL);
    start_time = timer.tv_sec + (timer.tv_usec/1000000.0); 

	// Create threads
	for (i = 0; i < num_threads; i++) {
		p = (int *) malloc(sizeof(int));
		*p = i;
		pthread_create(&threads[i], NULL, worker, (void *) (p));
	}

	// Join threads
	for (i = 0; i < num_threads; i++) {
		pthread_join(threads[i], NULL);
	}

	gettimeofday(&timer, NULL);
    end_time = timer.tv_sec + (timer.tv_usec/1000000.0); 

	// Compute the max_diff
	max_diff = max_diffs[0];
	for (i = 1; i < num_threads; i++) {
		if (max_diffs[i] > max_diff)
			max_diff = max_diffs[i];
	}

	printf("Maximum diff is: %f\n", max_diff);

    double tt = end_time - start_time;
	printf("Compute time for grid-size(%d), iterations(%d), "
            "cores(%d) is: %.3f seconds\n",
			grid_size, num_iters, num_threads, tt);

	free(grid);
	free(max_diffs);
	//free(arrive);

	exit(1);
}

